/*
 * Author: nicolas.bredeche(@)upmc.fr
 * Created on 2014-01-06 
 */

package picoevo.app.chemicalreactions;

import java.util.ArrayList;

import org.apache.commons.math3.exception.MaxCountExceededException;
import org.apache.commons.math3.exception.NoBracketingException;

import picoevo.core.Individual;
import picoevo.core.Population;
import picoevo.tools.Dictionary;
import picoevo.tools.Display;
import picoevo.tools.Utils;
import picoevo.tools.Version;

public class Individual_CR extends Individual {
	
	protected ArrayList<Reaction> _reactions = new ArrayList<Reaction>();

	protected double[][] matrix_Km = new double[Global_CR.enzymes.length][Global_CR.products.length];
	protected double[][] matrix_ki = new double[Global_CR.enzymes.length][Global_CR.products.length];

	public Individual_CR( Population __populationOwner) {
		super(__populationOwner);
		init();
	}
	
	public Individual_CR( Population __populationOwner, boolean __init)  // used by ODESimulation project (utils.IndividualGenerator)
	{
		super(__populationOwner);
		if(__init)
			init();
	}

	@Override
	public void init() 
	{
		// * initialize the K_m and k_i matrices for enzymatic reaction
		for ( int i = 0 ; i != Global_CR.enzymes.length ; i++)
		{
			for ( int j = 0 ; j != Global_CR.products.length ; j++)
			{
				matrix_Km[i][j] = getKmInitValue();
				matrix_ki[i][j] = getKiInitValue();  //Ftagn
			}
		}

		// * create reactions
		
		int nbReactions = (int) ( Math.random() * ( Global_CR.nbReactionsMax_atInit - Global_CR.nbReactionsMin_atInit ) + Global_CR.nbReactionsMin_atInit );

		for ( int i = 0 ; i != nbReactions ; i++ )
		{
			this.addReaction(new Reaction());
		}
		
		while ( this.isConsistent() == false ) // theoretically, an infinite loop is possible. In practical, it is unlikely to happen, ever.
		{
			this.forceConsistency();
			if ( _reactions.size() < Global_CR.nbReactionsMin_atInit )
			{ 
				this.addReaction(new Reaction());
			}
		} 
	}

	public void setMatrixKi ( double[][] matrix )
	{
		if ( matrix.length == matrix_ki.length && matrix[0].length == matrix_ki[0].length )
		{
			for ( int i = 0 ; i != matrix.length ; i++ )
				for ( int j = 0 ; j != matrix[0].length ; j++ )
					matrix_ki[i][j] = matrix[i][j];
		}
		else
			Display.critical("Individual_CR::setMatrixKi - dimensions do not match.");
	} //Ftagn
	
	public void setMatrixKm ( double[][] matrix )
	{
		if ( matrix.length == matrix_Km.length && matrix[0].length == matrix_Km[0].length )
		{
			for ( int i = 0 ; i != matrix.length ; i++ )
				for ( int j = 0 ; j != matrix[0].length ; j++ )
					matrix_Km[i][j] = matrix[i][j];
		}
		else
			Display.critical("Individual_CR::setMatrixKm - dimensions do not match.");
	}
	
	public double[][] getCopyOfMatrixKm ()
	{
		double[][] matrixCopy = new double[Global_CR.enzymes.length][Global_CR.products.length];
		for ( int i = 0 ; i != Global_CR.enzymes.length ; i++ )
			for ( int j = 0 ; j != Global_CR.products.length ; j++ )
				matrixCopy[i][j] = this.matrix_Km[i][j];
		return matrixCopy;
	}
	
	public double[][] getCopyOfMatrixKi ()
	{
		double[][] matrixCopy = new double[Global_CR.enzymes.length][Global_CR.products.length];
		for ( int i = 0 ; i != Global_CR.enzymes.length ; i++ )
			for ( int j = 0 ; j != Global_CR.products.length ; j++ )
				matrixCopy[i][j] = this.matrix_ki[i][j];
		return matrixCopy;
	}//Ftagn
	
	public void addReaction ( Reaction reaction )
	{
		this._reactions.add(reaction);
	}
	
	public void clearReactions()
	{
		this._reactions.clear();
	}
	
	private double getKmInitValue()
	{
		return Utils.getRandomPow(Global_CR.KmMinPowValue,Global_CR.KmMaxPowValue);
	} //Ftagn

	private double getKiInitValue()
	{
		return Utils.getRandomPow(Global_CR.KiMinPowValue,Global_CR.KiMaxPowValue);
	}
	
	// ---- ---- ---- ----
	
	private int checkConsistency()
	{
		int actions = 0;
		
		while ( this.isConsistent() == false || _reactions.size() < Global_CR.nbReactionsMin ) // theoretically, an infinite loop is possible. In practical, it is unlikely to happen, ever.
		{
			if ( _reactions.size() < Global_CR.nbReactionsMin )
			{ 
				this.addReaction(new Reaction());
				actions++;
			}
			actions += this.forceConsistency();
		} 
		
		return actions;
	}
	
	private boolean isConsistentWith( Reaction __reaction )
	{
		int[] signatureCandidate = __reaction.getSignature();
		
		for ( int i = 0 ; i != _reactions.size() ; i++ )
		{
			int[] signature = _reactions.get(i).getSignature();
			int score = 0;
			for ( int j = 0 ; j != signature.length ; j++ )
				if ( signatureCandidate[j] != signature[j] )
					score++;
			if ( score == signature.length) // ie. signature and signatureCandidate match perfectly
				return false;
		}
		
		return true;
	}

	private boolean isConsistent( )
	{
		int signatureLength = _reactions.get(0).getSignature().length;
		
		for ( int i = 0 ; i < _reactions.size()-1 ; i++ )
		{
			int[] signatureSource = _reactions.get(i).getSignature();

			for (int j = i+1 ; j < _reactions.size() ; j++ )
			{
				int[] signatureTarget = _reactions.get(j).getSignature();
				int score = 0;
				for ( int k = 0 ; k != signatureTarget.length ; k++ )
					if ( signatureSource[k] == signatureTarget[k] )
						score++;
				if ( score == signatureLength )
				{
					return false;
				}
			}
		}
		return true;
	}
	
	private int forceConsistency( )
	{
		int actions = 0;
		
		for ( int i = 0 ; i < _reactions.size()-1 ; i++ )
		{
			int[] signatureSource = _reactions.get(i).getSignature();

			for (int j = i+1 ; j < _reactions.size() ; j++ )
			{
				int[] signatureTarget = _reactions.get(j).getSignature();
				int score = 0;
				for ( int k = 0 ; k != signatureTarget.length ; k++ )
					if ( signatureSource[k] != signatureTarget[k] )
						score++;
				if ( score == 0 ) // check for duplicate
				{
					Display.warning("Individual_ChemReactions::forceConsistency - duplicate reaction is removed.");

					// remove one of the two duplicates
					if ( Math.random() < 0.5 )
					{
						_reactions.remove(i);
						continue;
					}
					else
						_reactions.remove(j);
					
					actions++;
				}
			}
		}
		
		return actions;
	}
	
	public Object clone ()
	{
		Individual_CR individualClone = null;
		
		individualClone = ( Individual_CR )super.clone();

		individualClone._reactions = new ArrayList<Reaction>();

		for ( int i = 0 ; i != this._reactions.size() ; i++ )
		{
			individualClone.addReaction( (Reaction)(this._reactions.get(i).clone()) );
		}
		
		individualClone.matrix_ki = new double[Global_CR.enzymes.length][Global_CR.products.length]; //Ftagn
		individualClone.matrix_Km = new double[Global_CR.enzymes.length][Global_CR.products.length];
		for ( int i = 0 ; i != Global_CR.enzymes.length ; i++ )
			for ( int j = 0 ; j != Global_CR.products.length ; j++ )
			{
				individualClone.matrix_ki[i][j] = this.matrix_ki[i][j]; //Ftagn
				individualClone.matrix_Km[i][j] = this.matrix_Km[i][j];
			}

		return individualClone;
	}




	public String getMatricesAsString()
	{
		String s_Km = "{ ";
		String s_ki = "{ ";//Ftagn
		for ( int i = 0 ; i != Global_CR.enzymes.length ; i++)
		{
			s_Km += "{ ";
			s_ki += "{ ";//Ftagn
			for ( int j = 0 ; j != Global_CR.products.length ; j++)
			{
				s_Km += ""+matrix_Km[i][j]+" ";
				s_ki += ""+matrix_ki[i][j]+" ";//Ftagn
			}
			s_Km += "}";
			s_ki += "}";//Ftagn
		}
		s_Km += " }";
		s_ki += " }";//Ftagn
		
		s_Km = "K_m["+Global_CR.enzymes.length+"]["+Global_CR.products.length+"] : " + s_Km + "\n";
		s_ki = "k_i["+Global_CR.enzymes.length+"]["+Global_CR.products.length+"] : " + s_ki + "\n";//Ftagn
		
		String s = s_Km + s_ki;//Ftagn
		
		return s;
	}
	
	@Override public void displayInformation() 
	{
		// general information
		super.displayInformation();
		
		// K_m and k_i values
		Display.info("# Constants wrt. enzymatic reaction:");
		Display.info_nocr(""+getMatricesAsString());
		
		// Reactions
		Display.info("# Reactions:");
		for ( int i = 0 ; i != _reactions.size() ; i++ )
			Display.info(_reactions.get(i).toString());
	}
	
	public String getDescription()
	{
		String s = "";
		s += "# individual \"" + this.getId() + "\"\n";
		
		// enzymes and products
		s += "enzymes:";
		for ( int i = 0 ; i != Global_CR.enzymes.length ; i++ )
			s += " " + Global_CR.enzymes[i];
		s += "\n";
		s += "products:";
		for ( int i = 0 ; i != Global_CR.products.length ; i++ )
			s += " " + Global_CR.products[i];
		s += "\n";
			
		
		// K_m and k_i values
		s += getMatricesAsString();
		
		// Reactions
		for ( int i = 0 ; i != _reactions.size() ; i++ )
			s += _reactions.get(i).toString() + "\n";
		
		return s;
	}

	private boolean mutateMatrixKm() 
	{
		this.matrix_Km[(int)Math.random()*this.matrix_Km.length][(int)Math.random()*this.matrix_Km[0].length] = getKmInitValue();
		return true;
	}

	private boolean mutateMatrixKi() 
	{
		this.matrix_ki[(int)Math.random()*this.matrix_ki.length][(int)Math.random()*this.matrix_ki[0].length] = getKiInitValue();
		return true;
	}//Ftagn

	private boolean mutateDecomposeReaction() // non-enzymatic reactions only 
	{	
		// count number of candidates
		int candidateCnt = 0;
		for ( int i = 0 ; i != this._reactions.size() ; i++ )
		{
			if ( this._reactions.get(i).leftEnzymeId == -1 )
				candidateCnt++;
		}
		
		// choose one 
		if ( candidateCnt > 0 )
		{
			int iSelected = -1;
			int iCnt = -1;
			int iCntTarget = (int)(Math.random()*candidateCnt);

			do {
				iSelected++;
				if ( this._reactions.get(iSelected).leftEnzymeId == -1 )
					iCnt++;
			} while ( iCnt != iCntTarget );
		
			// decompose 
			Reaction reaction = new Reaction();
			int[] signatureSource = this._reactions.get(iSelected).getSignature();
			
			// arbitrary initialize a new reaction until it is compatible (ie. consistent) with the source reaction) -- assume non-enzymatic reaction can be initialized
			// the last term of the while enables to avoid introducing null reaction.
			while ( 
					reaction.leftEnzymeId != -1 || 
					( reaction.rightProductOneId == signatureSource[1] && reaction.rightProductTwoId == signatureSource[2] ) || 
					( reaction.rightProductOneId == signatureSource[2] && reaction.rightProductTwoId == signatureSource[1] ) ||
					( reaction.rightProductOneId == signatureSource[4] && reaction.rightProductTwoId == signatureSource[5] ) ||
					( reaction.rightProductOneId == signatureSource[5] && reaction.rightProductTwoId == signatureSource[4] ) ||
					( reaction.rightProductOneId == -1 )
				  )  
			{
				reaction.init();
			}
			
			// update reactions
			int p1 = reaction.rightProductOneId;
			int p2 = reaction.rightProductTwoId;
			reaction.leftProductOneId = p1;
			reaction.leftProductTwoId = p2;
			reaction.rightProductOneId = this._reactions.get(iSelected).rightProductOneId; 
			reaction.rightProductTwoId = this._reactions.get(iSelected).rightProductTwoId;
			this._reactions.get(iSelected).rightProductOneId = p1;
			this._reactions.get(iSelected).rightProductTwoId = p2;
			
			// DEBUG
			/*
			int[] signature = this._reactions.get(iSelected).getSignature();
			//if ( signature[1] == -1 && signature[5] == -1 && signature[2] == -1 && signature[4] == -1 )
			{
				Display.debug("ERROR in decompose\nsignature:");
				for ( int i=0 ; i!=6 ; i++) 
					Display.debug("- signature["+i+"]: "+ signature[i]);
				Display.debug("");
				this.displayInformation();
			}
			/**/
			//DEBUG
			
			//FGJ
			/*
			* k+  k1+*k2+
			*-- = -------
			* k-  k1-*k2-
			*
			* mixingCoeff \in [0..1]
			*
			*k1x = sqrt(abs(kx)) * mixingCoeff
			*k2x = sqrt(abs(kx)) * (1 - mixingCoeff) 
			*
			*mixingCoeff = 0.5? au plus proche
                        * 
                        *La reaction est cependant plus lente?
                        * 
			*/
			
			double mixingCoeff = Math.random();
			double skp = Math.signum(this._reactions.get(iSelected).k_plus);
			double skm = Math.signum(this._reactions.get(iSelected).k_minus);

			double sqrtkp = Math.sqrt(this._reactions.get(iSelected).k_plus);
			double sqrtkm = Math.sqrt(this._reactions.get(iSelected).k_minus);

			reaction.k_plus = skp * sqrtkp * (1-mixingCoeff);
			reaction.k_minus= skm * sqrtkm * (1-mixingCoeff);

			this._reactions.get(iSelected).k_plus = skp * sqrtkp * mixingCoeff;
			this._reactions.get(iSelected).k_minus= skm * sqrtkm * mixingCoeff;

			if(skp < 0)
			{
				if(Math.random() >0.5)
				{
					this._reactions.get(iSelected).k_plus*=-1;
				}
				else
				{	
					reaction.k_plus *= -1;
				}
			}

			if(skm < 0)
			{
                if(Math.random() >0.5)
                {
                        this._reactions.get(iSelected).k_minus*=-1;
                }
                else
                {
                        reaction.k_minus *= -1;
                }

			}

			this.addReaction(reaction);
			
			return true;
			}

		Display.warning("Individual_ChemReactions::MutateReactionDecompose - no candidate to decompose");
		
		return false;
	}

	private boolean mutateMergeReaction() // non-enzymatic reactions only 
	{
		// count number of candidates
		int candidateCnt = 0;
		for ( int i = 0 ; i != this._reactions.size() ; i++ )
		{
			if ( this._reactions.get(i).leftEnzymeId == -1 )
				candidateCnt++;
		}
		
		// choose one 
		if ( candidateCnt >= 2 )
		{
			// search for non-enzymatic reaction
			// randomize array of indexes (to avoid ordering bias)
			// search for compatibility between pairs (candidate1.left = candidate2.right)
			// merge pair if found
			
			// create an array that points to possible candidates
			int[] candidateIndexes = new int[candidateCnt];
			for ( int i = 0, j = 0 ; i != this._reactions.size() ; i++ )
			{
				if ( this._reactions.get(i).leftEnzymeId == -1 )
				{
					candidateIndexes[j] = i; 
					j++;
				}
			}

			// randomize array (ie. enable to stop at the first couple of candidates, and still be sure that there is no bias)
			for ( int i = 0 ; i != candidateIndexes.length ; i++ )
			{
				int tmp = candidateIndexes[i];
				int j = (int) (Math.random()*candidateIndexes.length);
				candidateIndexes[i] = candidateIndexes[j];
				candidateIndexes[j] = tmp;
			}

			// scan for first possible pair of candidates
			
			//loop:
			for ( int i = 0 ; i != candidateIndexes.length-1 ; i++ )
			{
				for ( int j = i+1 ; j != candidateIndexes.length ; j++ )
				{
					int[] signature1 = this._reactions.get(candidateIndexes[i]).getSignature();
					int[] signature2 = this._reactions.get(candidateIndexes[j]).getSignature();
					if ( ( signature1[4] == signature2[1] && signature1[5] == signature2[2] ) ||  ( signature1[1] == signature2[4] && signature1[2] == signature2[5] ) )
					{
						if ( signature1[4] == signature2[1] && signature1[5] == signature2[2] )
						{
							//Display.debug("i: " + this._reactions.get(candidateIndexes[i]).toString());
							//Display.debug("j: " + this._reactions.get(candidateIndexes[j]).toString());
							signature1[4] = signature2[4]; 
							signature1[5] = signature2[5];
							if ( Reaction.checkConsistency(signature1) == true ) 
							{
								this._reactions.get(candidateIndexes[i]).setSignature(signature1);
								//Display.debug("R: " + this._reactions.get(candidateIndexes[i]).toString());
								this._reactions.remove(candidateIndexes[j]);
								return true;
							}
							//else
							//{
							//	Display.warning("Individual_ChemReactions::MutateReactionMerge - abort merge (new individual is empty)");
							//}
						}
						else
						{
							//Display.debug("j: " + this._reactions.get(candidateIndexes[j]).toString());
							//Display.debug("i: " + this._reactions.get(candidateIndexes[i]).toString());
							signature2[4] = signature1[4]; 
							signature2[5] = signature1[5];
							if ( Reaction.checkConsistency(signature2) == true ) 
							{
								this._reactions.get(candidateIndexes[j]).setSignature(signature2);
								//Display.debug("R: " + this._reactions.get(candidateIndexes[j]).toString());
								this._reactions.remove(candidateIndexes[i]);
								return true;
							}
							//else
							//{
							//	Display.warning("Individual_ChemReactions::MutateReactionMerge - abort merge (new individual is empty)");
							//	return false;
							//}
							/*
							 * 
							 * 
- Il y a surement plusieurs possibilité , mais tu peux par exemple prendre k+= min(k1+,k2+) et k-=k+*k1-*k2-/(k1+*k2+).
- Ca fais l’approximation que la transformation globale est limite par la reaction la plus lente… Dans tous les cas, tu aura du bruit ensuite pour ajuster tout ca...
- ne fait ca que si les reactions sont dans le bon sens (ie pas si tu as A<=>B et C<=>B), le sens étant donne par k+>k- (voir plus bas)


							 */
						}
						//break loop;
					}
				}
			}
		}
		else
			Display.warning("Individual_ChemReactions::MutateReactionMerge - no candidates to merge, or failed attempts (empty individual)");
		
		return false;
	}


	/*
	Mutation::individus
	- ajout d'une reaction (aleatoire, originale)
	- suppression d'une reaction (sous reserve que nb_reactions>nb_reactions_min)
	- decomposition d'une expression X=>Y en X=>Z et Z=>Y
	- recomposition d'une expression X=Y et Y=>Z en X=>Z

	Mutation::reaction
	- produit: ajout/suppression/mutation
	- enzyme: ajout/suppression/mutation (a droite et a gauche) ; Rq: ajout/suppression d'un enzyme a du sens p/r vitesse de reaction.
	
	returns: number of mutations (note that in some particular cases (failed attempt, or consistency forced), this is an upper bound)
	*/ 
	public int mutate() 
	{
		Display.info("Individual: " + this.getId());
	
		return mutate_onestepMutation();
		//return mutate_Continuous_parameters();
	}
	
	// ---- ---- ----
	
	/**
	 * Try to mutate only once. Note that this method may be called several times so as to travel for a given number of steps in the search space.
	 * @return number of actual mutation (which should be one)
	 */
	private int mutate_oneStepOneDimension()
	{
		
		
		return 0;
	}
	
	private int mutate_Continuous_parameters()
	{
		int mutationsIndiv = 0; //copier collé honteux
		int mutationsReactions = 0;


		int nbEnzymaticReactions = 0; 
		for ( int i = 0 ; i != this._reactions.size() ; i++ )
			if ( _reactions.get(i).isEnzymaticReaction() ) nbEnzymaticReactions++;
		
		double P_all;
		if ( nbEnzymaticReactions == 0 )
			P_all = Global_CR.P_onestep_indivParam;
		else
			P_all = Global_CR.P_onestep_indivParam + Global_CR.P_onestep_reactionParam;
		
		double dice = Math.random();


		if ( dice < ( Global_CR.P_onestep_indivParam / P_all ) ) // si aucune réaction ezimatique -> 1
		{
			/*int mutations = 0;
			do { 
				mutations = mutateIndividualParameters(); 
			} while ( mutations == 0 );
			mutationsIndiv += mutations;*/
			int dice2 = (int)((Math.random())*(this._reactions.size() - nbEnzymaticReactions));
			int j = -1;
			for ( int i = 0 ; j != dice2 ; i++ )
				if ( !this._reactions.get(i).isEnzymaticReaction() ) j++;
			this._reactions.get(j).mutateParameters();

			mutationsReactions++;
		}
		else 
		{
			while( mutateIndividualParameters() == 0) ;

			// mutate one reaction parameters (should not happen if no enzymatic reactions exist.
			/*int dice2 = (int)( Math.random() * nbEnzymaticReactions );
			int j = -1;
			for ( int i = 0 ; j != dice2 ; i++ )
				if ( this._reactions.get(i).isEnzymaticReaction() ) j++;
			this._reactions.get(j).mutateParameters(); //mutate parameters => k+ k- => non enzymatic aussi? alors que mutateIndividualParameters km ki...*/

			mutationsIndiv++;
		}
		
		
		
		
		return (mutationsIndiv+mutationsReactions);
	}
	
	private int mutate_onestepMutation()
	{
		// * individual-level mutation
		
		int mutationsIndiv = 0;
		int mutationsReactions = 0;

		// check for enzymatic reactions (impact which mutations are possible)
		
		int nbEnzymaticReactions = 0; 
		for ( int i = 0 ; i != this._reactions.size() ; i++ )
			if ( _reactions.get(i).isEnzymaticReaction() ) nbEnzymaticReactions++;
		
		double P_all;
		if ( nbEnzymaticReactions == 0 )
			P_all = Global_CR.P_onestep_indivParam + Global_CR.P_onestep_indivStructure + Global_CR.P_onestep_reactionStructure;
		else
			P_all = Global_CR.P_onestep_indivParam + Global_CR.P_onestep_indivStructure + Global_CR.P_onestep_reactionParam + Global_CR.P_onestep_reactionStructure;

		// pick one mutation type

		double dice = Math.random();

		if ( dice < ( Global_CR.P_onestep_indivParam / P_all ) )
		{
			int dice2 = (int)((Math.random())*(this._reactions.size() - nbEnzymaticReactions));
			int j = -1;
			for ( int i = 0 ; j != dice2 ; i++ )
				if ( !this._reactions.get(i).isEnzymaticReaction() ) j++;
			this._reactions.get(j).mutateParameters();

			mutationsReactions++;
			/*int mutations = 0;
			do { 
				//mutations = mutateIndividualParameters(); //Uniquement pour les enzy !
				
			} while ( mutations == 0 );
			mutationsIndiv += mutations; */
		}
		else
			if ( dice < ( Global_CR.P_onestep_indivParam + Global_CR.P_onestep_indivStructure ) / P_all )
			{
				int mutations = 0;
				do { 
					mutations = mutateIndividualStructure(); 
				} while ( mutations == 0 );
				mutationsIndiv += mutations; 
			}
			else
				if ( dice < ( Global_CR.P_onestep_indivParam + Global_CR.P_onestep_indivStructure + Global_CR.P_onestep_reactionStructure ) / P_all )
				{
					int dice2 = (int)( Math.random() * this._reactions.size() );
					int mutate = 0, cnt = 0;
					do {
						mutate = this._reactions.get(dice2).mutateStructure(this);
						cnt++;
						if ( cnt == 1000 )
							Display.error("Individual_CR::mutate_onestepMutation() -- failed to mutate after " + cnt + " tentatives. Aborting."); // It is very unexpected that this should ever occur.
					} while ( mutate == 0 ); // force at least one mutation -- theoretically an infinite loop is possible

					mutationsReactions++;
				}
				else
				{
					mutationsIndiv += mutateIndividualParameters();
					// mutate one reaction parameters for non-enzymatic reactions (should not happen if no enzymatic reactions exist).
					//int dice2 = (int)( this._reactions.size() - Math.random() * nbEnzymaticReactions ); //ftagn FGJ Et tout le reste !
//					int dice2 = (int)(Math.random())*(this._reactions.size() - nbEnzymaticReactions);
//					int j = -1;
//					for ( int i = 0 ; j != dice2 ; i++ )
//						if ( !this._reactions.get(i).isEnzymaticReaction() ) j++;
//					this._reactions.get(j).mutateParameters();
//
//					mutationsReactions++;
				}
					
		// * ensure consistency
		
		checkConsistency();
		
		// * display information
		
		Display.info("Individual-level: " + mutationsIndiv + " mutation(s).");
		Display.info("Reaction-level:   " + mutationsReactions + " mutation(s).");
		
		if ( mutationsIndiv + mutationsReactions == 0 )
			Display.warning("Individual_CR::mutate_onestepMutation() -- Mutation failed.");
		
		return ( mutationsIndiv + mutationsReactions );
	}
	
	// ---- ---- ----
	
	/**
	 * Try to mutate everything.
	 * @return number of actual mutations
	 */
	private int mutate_defaultMethod()
	{
		// * individual-level mutation
		
		int mutationsIndiv = 0;

		mutationsIndiv += mutateIndividualParameters();
		mutationsIndiv += mutateIndividualStructure();
		
		// * reaction-level mutation

		int mutationsReactions = 0;
		int cnt = 0;
		
		do {
			for ( int i = 0 ; i != this._reactions.size() ; i++ )
			{
				Reaction reaction = this._reactions.get(i);
				mutationsReactions += reaction.mutateStructure(this);
				if ( reaction.isEnzymaticReaction() )
					mutationsReactions += reaction.mutateParameters();
			}
			
			cnt++;
			if ( cnt == 1000 )
			{
				Display.error("Individual_CR::mutate_onestepMutation() -- failed to mutate after " + cnt + " tentatives. Aborting."); // It is very unexpected that this should ever occur.
			}
		} while ( mutationsReactions == 0 ); // force at least one mutation -- theoretically an infinite loop is possible
		
		// * ensure consistency
		
		checkConsistency();
		
		// * display information
		
		Display.info("Individual-level: " + mutationsIndiv + " mutation(s).");
		Display.info("Reaction-level:   " + mutationsReactions + " mutation(s).");
		
		return mutationsIndiv + mutationsReactions;
	}
	
	
	// ---- ---- ----
	
	private int mutateIndividualParameters() // => one mutation at most.
	{
		int mutationsIndiv = 0;
		
		// mutate values in matrices
		
		double dice = Math.random();
		double P_all = Global_CR.P_mutate_matrix_Km + Global_CR.P_mutate_matrix_ki;  
		
		if ( dice < Global_CR.P_mutate_matrix_Km / P_all )
		{
			Display.info("Individual: mutate K_m.");
			mutationsIndiv++;
			mutateMatrixKm();
		}
		else		
		{
			Display.info("Individual: mutate K_i.");
			mutationsIndiv++;
			mutateMatrixKi(); //Ftagn
		}
		
		return mutationsIndiv;
	}

	
	// ---- ---- ----

	private int mutateIndividualStructure() // => one mutation at most.
	{
		int mutationsIndiv = 0;
		
		double dice = Math.random();
		double P_all = Global_CR.P_mutate_reaction_add + Global_CR.P_mutate_reaction_remove + Global_CR.P_mutate_reaction_decompose + Global_CR.P_mutate_reaction_merge;
		double P_current = 0;
		
		// add a random reaction
		P_current += Global_CR.P_mutate_reaction_add;
		if ( dice < P_current / P_all && this._reactions.size() < Global_CR.nbReactionsMax )
		{
			Display.info("Individual: add a reaction.");
			mutationsIndiv++;
			this.addReaction(new Reaction());
		}
		else
		{
			// remove a reaction
			P_current += Global_CR.P_mutate_reaction_remove;
			if ( dice < P_current / P_all && this._reactions.size() > Global_CR.nbReactionsMin )
			{
				Display.info("Individual: remove a reaction.");
				mutationsIndiv++;
				_reactions.remove( (int)(Math.random()*_reactions.size()-1) );
			}
			else
			{
				// decompose a reaction (non-enzymatic only)
				P_current += Global_CR.P_mutate_reaction_decompose;
				if ( dice < P_current / P_all && this._reactions.size() < Global_CR.nbReactionsMax )
				{
					Display.info("Individual: decompose a reaction.");
					mutationsIndiv++;
					mutateDecomposeReaction();
				}
				else
				{
					// merge two compatible reactions (non-enzymatic only)
					//P_current += Global_CR.P_mutate_reaction_merge; // useless as it sums to 1.0
					if ( this._reactions.size() > Global_CR.nbReactionsMin )
					{
						Display.info("Individual: merge two reactions.");
						mutationsIndiv++;
						mutateMergeReaction();
					}
				}
			}
		}
	
		return mutationsIndiv;
	}
	
	// ---- ---- ----
	
	@Override
	public void evaluate() // TODO: plug into chemical reactions simulator
	{
		double fitnessValue = -1;
		
		//String s = this.getDescription(); // s now contains the comprehensive description of this individual, ready to be evaluated by Nath.
		
		//fitnessValue = Math.random(); // pure random
		//fitnessValue = this._populationOwner.getGenerationNumber(); // monotonous increase
		
		/**/
        double startTime = System.currentTimeMillis();
		try 
		{
			//fitnessValue = odesimulation.OdeImplicitSolver_Nicolas.getFitnessValue_1( this );
			fitnessValue = odesimulation.OdeImplicitSolver_Nicolas.getFitnessValue_2( this.getDescription() );			
		}
		catch ( MaxCountExceededException e )
		{
			if ( this._populationOwner.getProperties().getStringProperty("OptimisationFlag").compareTo(Boolean.toString(Dictionary.Minimization)) == 0 )
				fitnessValue = Double.MAX_VALUE;
			else
				fitnessValue = -Double.MAX_VALUE;
			Display.info("Fitness value: " + fitnessValue + " (discarded evaluation)");
		}
		catch ( NoBracketingException e )
		{
			if ( this._populationOwner.getProperties().getStringProperty("OptimisationFlag").compareTo(Boolean.toString(Dictionary.Minimization)) == 0 )
				fitnessValue = Double.MAX_VALUE;
			else
				fitnessValue = -Double.MAX_VALUE;
			Display.info("Fitness value: " + fitnessValue + " (discarded evaluation)");
		}
        Display.info("(chrono: "+ ((System.currentTimeMillis()-startTime)/1000) +"s).\n");
		/**/
		
		this.setFitness(fitnessValue); 
	}

}

